#!/usr/bin/python
'''
Translates SAM to HAM format.
'''
import sys
import ham
import h5py
import numpy as np
import time
import logging

# parameters.
sam_file = sys.argv[1]
ham_file = sys.argv[2]
#num_threads = int(sys.argv[3])

#buffer = 32212254720		# 30GB
#buffer = 10737418240		# 10GB
buffer = 1073741824		# 1GB
#buffer = buffer / 100
logging.basicConfig(level=logging.DEBUG, format='[%(levelname)s] (%(threadName)-10s) %(message)s', )



# setup simple structures.
logging.debug("Setting up data structures.")
ham.ham_dt = ham.ham_50_dt
ham_obj = ham.HamFile(ham_file, cprs=True)
lines = []
sz = 0
i = 0
j = 0
chunk_cnt = 0
init_time = time.time()


# read in sam file in chunks.
fin = open(sam_file, "rb")
while 1 == 1:

	# read in new chunk.
	logging.debug("reading chunk %i" % chunk_cnt)
	lines = fin.readlines(buffer)
	chunk_cnt += 1
		
	
	# check if we are done.
	if lines == []: break
	
	# instantiate a ham array.
	sz = len(lines)
	logging.debug("creating ham chunk of size %i." % sz)
	ham_data = ham.create_array(sz)
	
	# loop over chunk.
	logging.debug("processing chunk.")
	for i in range(sz):
		
		# skip headers.
		if lines[i][0] == "@": continue
				
		# process required part of line.
		tmp = lines[i].strip().split("\t")

		ham_data[i]['QNAME'] = tmp[0]
		ham_data[i]['FLAG'] = int(tmp[1])
		ham_data[i]['RNAME'] = tmp[2]
		ham_data[i]['POS'] = int(tmp[3])
		ham_data[i]['MAPQ'] = int(tmp[4])
		ham_data[i]['CIGAR'] = tmp[5]
		ham_data[i]['MRNM'] = tmp[6]
		ham_data[i]['MPOS'] = int(tmp[7])
		ham_data[i]['ISIZE'] = int(tmp[8])
		ham_data[i]['SEQ'] = tmp[9]
		ham_data[i]['QUAL'] = tmp[10]
	
		
		# throw the rest into optional.
		if len(tmp) > 10:
			ham_data[i]['OPTIONAL'] = '\t'.join(tmp[11::])

	# delete chunk.
	del lines		
	
	# save chunk of ham obj.
	ham_obj.save_array(ham_data)

# close input file.
fin.close()

# close ham pointer.
ham_obj.close()
	
# done.
logging.debug("total time: %f" % (time.time() - init_time))
